Metabarcoding of canopy arthropods reveals negative impacts of forestry insecticides on community structure across multiple taxa

Insect outbreaks pose a major challenge in the management and development of land resources. In cropping systems, these eruptions are frequent as relatively low diversity transient nature agroecosystems hinders trophic regulation herbivores population (Letourneau, 2012). contrast, insect forests characterised mainly by sporadic often spatially synchronous eruptions, resulting from combination abiotic biotic processes whose relative importance is still poorly understood 2012; Liebhold et al., Although infrequent, forest insects can substantially impede timber production (MacLean, 2016). Lepidopteran folivores, such gypsy moth Lymantria dispar L. (Lepidoptera: Erebidae), notably responsible for growth reduction increased tree mortality temperate (Lobinger, 1999; MacLean, 2016; Twery, 1991) routinely controlled aerial applications insecticides. By contrast with systematic preventative approach that predominates treatments applied curatively to suppress outbreaking populations before defoliation becomes significant (Thompson, 2011). However, despite occasional operations selectivity approved substances, use insecticides has raised concerns about its side effects on arthropod fauna, particularly non-target Lepidoptera (Schweitzer, 2004; Severns, 2002). Central European deciduous oak woodlands, insecticide generally late-April mid-May against spring feeding defoliators (Gößwein & Lobinger, 2014). During this period, larvae make up most biomass crowns cornerstone food web (Southwood 2004). suppressing populations, alter species interactions competition (Leroy, Gossner, 2021) predation (Sample 1993), potentially triggering complex chain feedback loops across entire communities (Fleeger, 2020). while direct impacts forestry well-understood, indirect more difficult study field remain known. Studies aerially sprayed suffer several methodological limitations. First, sampling techniques biased towards certain taxa due life history habitat among species, making it sample comprehensively (Ozanne, 2005). Second, taxonomic expertise rarely available all sampled (Brown, 2005), assessments usually focused few well-known rather than diverse communities. Last, primarily tested coarse ranks (often order or family) groupings reflecting hypotheses sensitivity species. Most products defoliator non-systemic ingestion moulting inhibitor diflubenzuron (Marx, 1977) Lepidoptera-specific gut poison Bacillus thuringiensis var. kurstaki (Gill 1992). Thus, arthropods be classified into hypothetical groups based their habits. For example, chewing expected sensitive leaf miners sucking feed internal plant tissues. The analyses may reveal toxic substances which would masked when ecologically investigated. This help improve our understanding forests. Recent advances amplicon-based metabarcoding have brought new opportunities circumvent shortcomings determination morphological examination (Hebert Gregory, With implementation next generation sequencing (NGS) within framework, millions amplified cytochrome c oxidase 1 (CO1) genes now sequenced simultaneously searched reference barcode libraries bulk-identify presence raw samples (Baird Hajibabaei, last decade, rapid sizeable collaborative databases Barcode Life Data System (BOLD; Ratnasingham Hebert, 2007) greatly scope approaches biodiversity. BOLD, similar CO1 DNA barcodes clustered attributed unique identifier (Barcode Index Number, BIN) closely resemble used proxy resolution (Ratnasingham 2013). Delineations genetic distance unbiased asymmetric state expertise, enabling comparison studies. approach, combined comprehensive canopy pyrethrum knockdown offers robust framework assess patterns community structure not devoid flaws. inferring abundances sequence read pivotal issue current approaches. Sequence divergence differences copy number targeted locus, well variation body mass specimen, lead amplification hindering detection some estimation abundance reads largely inaccurate (Elbrecht 2017; Krehenwinkel 2017). present study, we assessed short-term (5 weeks post-application) DFB BTK control early summer woodlands. We maximise implemented simple protocol integrating order-level sorting subsampling obtain reliable estimates numbers. then coarsely ecological illustrating exposure knowledge toxicological environmental dynamics species' life-history traits. Specifically, distinguished affected via primary (externally [DFB; BTK], other externally [DFB only]), secondary (parasitoids, predators, scavengers only]) unaffected (internally scavengers). trees treated DFB, left unsprayed, test whether reliably predicted mere basis effect hypotheses. experiment was conducted three neighbouring oak-dominated stands, ‘Bauernschlag’ (52 ha; 50.05706°N, 10.08311°E), ‘Vorberg’ (58 50.0534°N, 10.08899°E) ‘Brunnholz’ (6 50.03170°N, 10.07715°E), region Schweinfurt (Lower Franconia, Bavaria, Germany) 2017. sites were selected surveys during previous winter (October 2016–February 2017) predicting patchy densities ranging endemic potential outbreak levels. Within 10 structurally homogeneous areas stands (blocks), delimited two 3,000-m2 (100 × 30 m) treatment plots randomly assigned one (DFB BTK) plot. Plots separated buffer zones at least 40-m wide, considered safe prevent off-target spraying air conditions (German Federal Office Consumer Protection Food Safety, 2019). A detailed description experimental design provided Appendix S1. commercial formulations licensed DiPel® ES (BTK strain ABTS 351 [HD-1]; 33.2 g/L; 17,600 IU/mg; Cheminova Deutschland) Dimilin® 80 WG (800 g DFB/kg; Spiess-Urania Chemicals) used. Both selective They non-systemic, means they surfaces plants after application, do translocate tissues must ingested effective (Bull Ivie, 1978; Gill 1992; Grosscurt, 1978). Their main difference lies spectrum action, being specific Lepidoptera, juvenile stages upon oral contact. Table summarises ecotoxicological data toxicity both substances. maximal legal rate 3 L/ha 75 g/ha respectively. administered Bell 47 helicopter 17 May 2017 between 07:00 09:00 dry conditions. accuracy coverage water-sensitive spray cards (agrotop, Obertraubling, placed floor intervals 2–3 m along straight lines connecting centroid adjacent plots. total 131 exposed transect crossing four blocks (Appendix S1) collected replaced each prey flight inspect droplet deposition patterns. Uniform observed plot boundaries. Off-target spraying, measured droplets zones, did occur beyond away (in case fogging) free-living crowns. One mature (i.e. DBH > 35 cm) close centre Patchy egg masses variable levels blocks: (blocks 1, 2 10), moderate 4–9) severe (block 11; S1). hence area block compare high To measure efficacy treatments, 10-cm-wide band barricade tape around trunk focal tree, allowed us late-instar taking shelter daytime. caterpillars counted 12 13 June single occasion 19 22 2017, 5 corresponds peak (and defoliation). Swingfog® SN-50 fogging machines (Swingtec GmbH) 1% suspension petroleum white oil generate thermal cloud carries canopy. Fogging dusk dawn efficiency 2005) simultanesouly block) nullify impact drift fog machine operated until crown coated fog, determined visually using electric torch lights. time ranged min per tree. fallen 15-m2 tarpaulin sheets beneath stored them 99% ethanol further analysis (one sheet tree). Permits granted Regional Council Lower Franconia (Regierungspräsidium Unterfranken) §45 para 7 cl no Nature Conservation Act. No ethical approval required study. Adults heterometabolous undergoing incomplete metamorphosis adults juveniles morphologically similar) sorted 28 (mostly suborder, referred ‘sorted taxa’; see Figure S1 complete list). Holometabolous radically differ appearance) groups: Lepidoptera/Symphyta prolegs present) Coleoptera/Neuroptera/Raphidioptera absent). Because already stems days canopy, individuals excluded analyses. Based information relevant systematics taxa, taxon ‘optimal’ method: (a) expert taxonomist (MPH) (b) through NGS. MPH group included Araneae, Opiliones, Pseudoscorpiones, adult Coleoptera, Orthoptera, Myriapoda Isopoda. vast majority oak-associated described easily identifiable experts network. deemed challenging NGS group, owing various obstacles, many cryptic (e.g. Hemiptera), dependence non-ubiquitous characters male genitalia families Diptera Neuroptera), Hymenoptera) lack Psocodea). derive assemblage barcoded, first differentiated (Figure divided equal-size fractions, subsamples taxon, individual distributed numbered vials, starting random henceforth proceeding sequentially, individual, ascending 1). necessary because always multiple vial ensured there bias what receive higher individuals. aimed increase capacity, reducing amount sample. It also an aggregating incidence level sum presence/absence subsamples), non-null (only subsample) (present subsamples). Each 480 filled preserve DNA. standard methods extraction, amplification, NGS, assignment BIN putative species; BINs include undescribed subsequent full pipeline (DNA bioinformatic processing) detail S2. filtered out bias-corrected below maximum subsamples, likely contaminants. evaluate performance calculated proportion observable mismatches (false false absence taxa) cross-checking matrix obtained counting metabarcoding; Section 2.4.1) pooling taxa. compared filtering. estimate converted presence–absence aggregated data. value 0 method, fitted linear regressions strength relationship tree) proxies step protocol, namely contaminant filtering), correction filtering) (c) latter, models influence estimation: (4 subsample (12 scoring physiological sensitivity. single-species system, defined mode action insecticides: physiologically whereas only regarded BTK. Ecological possibility ecology. As truly exposure, surface tissue epiphyte grazers) directly concealed shelters herbivores, miners, shoot wood borers) avoid exposure. Predators, parasitoids exposed, come contact residues contaminated host. Importantly, non-predatory known association any overstorey understorey recorded ‘tourists’ sense Moran Southwood, 1982) analyses, very sporadic. Fully borers stage but been active acute phase following application). Hence, decided analysis, regardless. Saproxylic dead associated fungi) wood-boring beetles live trees) grouped herbivore guild plasticity identical predictions regarding scores sensitivities, drafted six so-called (tourists included; 2). Feeding trait retrieved mining specialised literature (see S3 bibliography). traits highly conserved, family systematically dropping matched seems assume braconids system parasitoids, regardless names). biologically families, phytophagous unassigned dropped. contribution different orders shown statistical procedures separately identified (further assemblage) type, BINs. clarity, will refer metric type. average common (?-diversity), mean (?-diversity) compositional heterogeneity (?-diversity) whole five lepidopteran external feeders, non-lepidopteran predators). Only (represented exclusively xylophagous beetles) predators sufficiently represented analysed community. dosing mistake night knockdown, (16 too under-sampled Therefore, remaining (24 trees). generalised mixed-effect model (both dummy-coded) fixed intercept nested site. same investigate quantitative response treatments. shortcoming method changes atmospheric affect cloud, leading considerable trees. took note large extent computing rarefaction curves S3). Abundance-based allows correcting uneven effort standardising smallest size (Gotelli Colwell, 2001). driven true muted process. driven, instead standardised coverage, asymptotic (Chao Jost, extrapolated richness twice abundance-based extrapolation curves. lowest ‘base coverage’ Chao Jost (2012). Species interpolated corresponding base size’). Following estimated (?-diversity) sizes treatment. treatment, halved correct doubled replicates (n = n 6 BTK). identify responses looked species–treatment associations performing indicator (De Cáceres 2010) involving 67 ?2 treatment). point-biserial correlation coefficient (rpb) corrected unequal sizes. significance (P) two-sided permutations tests (1,000 permutations), p-values subsequently adjusted testing Šidák correction. All performed R 4.0.3 (R Core Team, Mixed likelihood estimations ‘glmmTMB’ (Brooks diagnosed package ‘DHARMa’ (Hartig, (p < 0.05) type II Wald function Anova ‘car’ (Fox Weisberg, Interpolation ?- ?-diversity ?-diversity multiplicative partitioning ‘iNEXT’ (Hsieh uncertainty constructing 95% bootstrap confidence (500 replicates) ‘iNEXT’. Confidence (family-wise error 18 FWER 0.28%). statistically 5% indicated non-overlapping ‘indicspecies’ Legendre, 2009). dataset 13,634,135 9,042 OTUs 2,501 85 matches COI sequences <500 bp ambiguous bases, pending analyses). reduced 2,334 BIN-species discordances S2). total, 3,667 (11,563,819 reads) 1,214 retained removal non-sorted human contamination), pairwise identity 97% (<0.005% sample). After numbers filtering (<1% BIN), final consisted 10,846,068 reads. Before reads, output 5.6% positives contaminants artefacts) 4.7% negatives undetected diversity; 2, a1) Small-sized Acari (36.8% negatives), Thysanoptera (12.2%) Sternorrhyncha (11.8%) had failures detections larger sized fraction down 0.8% 8.3% a2). False slightly (+12.5%), Collembola (+8.3%) grouping Coleoptera larvae, Neuroptera Raphidioptera (+6.6%; Metabarcoding 25,949 processed delineation metabarcoding. 1,213 detected taxonomy level, 1,145 (94%) 960 (79%) genus 757 (62%) level. Among percentage varied 96% (113 BINs) 22% (23 BINs). species-level 48 Hymenoptera 63 74% Hemiptera respectively b1). Morphology 4,171 4,771 successfully These (Coleoptera; 2,715 2,717), orthopterans (275), harvestmen (Opiliones; 14), pseudoscorpions (14), myriapods (36) isopods (3), 65% spiders (Araneae; 1,114 1,712). comprised almost could (27; 2%) (571; 34%; b2). weak (t[22] 0.93, p 0.364; R2 0.038; c1). Scaling marginally 1.46, 0.160; 0.088; c2). 11.08, 0.001; 0.810; c3) 13.56, 0.860; c4) strongly linearly related abundance. Increasing decreased shrinkage estimates, slope 6.74 units) 4.79 units). under bands (5- 43-fold respectively; ?2(2) 60.01, Gypsy 13% although pest effective. (Table Insecticide overall 4). lower controls groups. pronounced (55% fewer species) (25% DFB. trend group. pattern loss approximately (?-diversity; 4), suggesting 52% 27% DFB-treated re